home *** CD-ROM | disk | FTP | other *** search
/ Space & Astronomy / Space and Astronomy (October 1993).iso / mac / programs / space / AA_51.ZIP / OMERCURY.C < prev    next >
C/C++ Source or Header  |  1993-02-13  |  2KB  |  90 lines
  1. /* Orbital elements and perturbations for the planet Mercury
  2.  * using formulas given by Meeus.
  3.  * All values computed are for equinox of date.
  4.  */
  5.  
  6.  
  7.  
  8. #include "planet.h"
  9. #include "kep.h"
  10.  
  11. /* Calculate orbital elements for Mercury
  12.  */
  13. int omercury(e,J)
  14. struct orbit *e;
  15. double J;
  16. {
  17. /* Orbital elements referred to equinox of date: */
  18. e->epoch = J;
  19. /* calculate all the mean anomalies */
  20. manoms(J);
  21. e->M = M1;
  22.  
  23. /* semimajor axis */
  24. e->a = 0.3870986;
  25. /* eccentricity */
  26. e->ecc = ( -0.000000030*T + 0.00002046)*T + 0.20561421;
  27.  
  28. #if OFDATE
  29. e->equinox = J;
  30. /* inclination */
  31. e->i = ( -0.0000183*T + 0.0018608)*T + 7.002881;
  32. /* argument of the perihelion */
  33. f = ( 0.0001208*T + 0.3702806)*T + 28.753753;
  34. e->w = mod360(f);
  35. /* longitude of the ascending node */
  36. f = ( 0.0001739*T + 1.1852083)*T + 47.145944;
  37. e->W = mod360(f);
  38. /* mean longitude */
  39. f = ( 0.0003011*T + 149474.07078)*T + 178.179078;
  40. e->L = mod360(f);
  41. #else
  42. e->equinox = J2000;
  43. /* inclination */
  44. e->i = ((-3.5e-8*T + 6.9e-7)*T - 0.0059556)*T + 7.010678;
  45. /* argument of the perihelion */
  46. f = ((4.3e-8*T + 7.445e-5)*T + 0.2842765)*T + 28.839814;
  47. e->w = mod360(f);
  48. /* longitude of the ascending node */
  49. f = ((-6.8e-8*T - 8.844e-5)*T - 0.1254715)*T + 48.456876;
  50. e->W = mod360(f);
  51. /* mean longitude */
  52. f = e->M + e->w + e->W;
  53. e->L = mod360(f);
  54. #endif
  55. return(0);
  56. }
  57.  
  58.         
  59. /* corrections due to perturbations
  60.  * These are added in after solving Kepler's equation.
  61.  */
  62. int cmercury(e)
  63. struct orbit *e;
  64. {
  65. double p;
  66.  
  67. /* Mercury's perturbations in longitude: */
  68. f = (5.0*M2 - 2.0*M1 + 12.220)*DTR;
  69. p = 0.00204*cos(f);
  70. f = (2.0*M2 - M1 - 160.692)*DTR;
  71. p += 0.00103*cos(f);
  72. f = (2.0*M5 - M1 - 37.003)*DTR;
  73. p += 0.00091*cos(f);
  74. f = (5.0*M2 - 3*M1 + 10.137)*DTR;
  75. p += 0.00078*cos(f);
  76. e->L += p*DTR;
  77.  
  78. /* Mercury's perturbations in radius vector */
  79. f = (2.0*M5 - M1 + 53.013)*DTR;
  80. p = 0.000007525*cos(f);
  81. f = (5.0*M2 - 3.0*M1 - 259.918)*DTR;
  82. p += 0.000006802*cos(f);
  83. f = (2.0*M2 - 2.0*M1 - 71.188)*DTR;
  84. p += 0.000005457*cos(f);
  85. f = (5.0*M2 - M1 - 77.75)*DTR;
  86. p += 0.000003569*cos(f);
  87. e->r += p;
  88. return(0);
  89. }
  90.